[数据分析与可视化] Python绘制数据地图1 |
您所在的位置:网站首页 › python 开源数据分析 › [数据分析与可视化] Python绘制数据地图1 |
本文主要介绍GeoPandas的基本使用方法,以绘制简单的地图。GeoPandas是一个Python开源项目,旨在提供丰富而简单的地理空间数据处理接口。GeoPandas扩展了Pandas的数据类型,并使用matplotlib进行绘图。GeoPandas官方仓库地址为:GeoPandas。GeoPandas的官方文档地址为:GeoPandas-doc。 本文所有代码见:Python-Study-Notes GeoPandas推荐使用Python3.7版本及以上,运行环境最好是linux系统。GeoPandas安装命令如下: pip install geopandas 如果上述命令安装出问题,则推荐使用conda安装GeoPandas,命令如下: conda install geopandas 或: conda install --channel conda-forge geopandas # jupyter notebook环境去除warning import warnings warnings.filterwarnings("ignore") # 查看geopandas版本 import geopandas as gpd gpd.__version__ '0.10.2' 目录1 基础结构1.1 GeoSeries1.1.1 GeoSeries基础元素1.1.2 GeoSeries属性与方法1.2 GeoDataFrame2 数据集基本操作2.1 数据集读写2.2 数据集筛选2.2.1 Bounding Box Filter2.2.2 Geometry Filter2.2.3 Row Filter2.2.4 cx索引器3 坐标参考系3.1 坐标参考系简介3.2 EPSG编码3.3 在GeoPandas中管理坐标参考系3.3.1 坐标参考系管理3.3.2 选择合适的投影坐标参考系4 参考4.1 GeoPandas使用4.2 相关资料 1 基础结构GeoPandas在Pandas的基础上增加对地理空间数据的支持。其主要结构为GeoSeries和GeoDataFrame。如下图所示,整张图片表示一个GeoDataFrame,行表示地图上的某一区域,列表示该区域的属性数据,其中geometry表示该区域的位置信息。每一个区域的位置形状用一个或多个几何对象表示,几何对象在GeoPandas中以GeoSeries格式数据构成。 在Geogeopandas中,GeoSeries用于构成基础的几何对象。GeoSeries由Python库Shapely实现的几何对象构成,这些几何元素可以是点(集)、线(集)、面(集): Points / Multi-Points Lines / Multi-Lines Polygons / Multi-Polygons 在jupyter notebook中, 可以直接以svg格式展示GeoSeries中的单个元素,使用如下: Points 点 # Point(x,y)用于创建单个点 from shapely.geometry import Point p = gpd.GeoSeries([Point(1, 1), Point(2, 2), Point(3, 3)],index=['p1','p2','p3']) p p1 POINT (1.00000 1.00000) p2 POINT (2.00000 2.00000) p3 POINT (3.00000 3.00000) dtype: geometry # 展示单个元素 p[0]Multi-Points 点集 # MultiPoint([(x1,y1), ... ,(xn,yn)])用于创建多个点构成的集合 from shapely.geometry import MultiPoint mp = gpd.GeoSeries([MultiPoint([(0, 1), (1, 0)]), MultiPoint([(0, 0), (1, 1)])], index=['mp1', 'mp2']) mp mp1 MULTIPOINT (0.00000 1.00000, 1.00000 0.00000) mp2 MULTIPOINT (0.00000 0.00000, 1.00000 1.00000) dtype: geometry # 展示单个元素 mp[0]Lines 线 # LineString([(x1,y1), ... , (xn,yn)])用于创建单条线段 from shapely.geometry import LineString l = gpd.GeoSeries([LineString([(0.5, 0), (1, 2), (0, 1)]), LineString([(0, 0), (3.3, 1)])], index=['l1', 'l2']) l l1 LINESTRING (0.50000 0.00000, 1.00000 2.00000, ... l2 LINESTRING (0.00000 0.00000, 3.30000 1.00000) dtype: geometry # 展示单个元素 l[0]MultiLineString 线集 # MultiLineString([LineString1,LineString2])用于创建多条线段的集合 from shapely.geometry import MultiLineString ml = gpd.GeoSeries([MultiLineString([[(1, 0), (2, 1)], [(-1, 0), (0, 5), (1, 0)]]), MultiLineString([[(0.5, 0), (1, 2), (0, 1)]])], index=['m1','m2']) ml m1 MULTILINESTRING ((1.00000 0.00000, 2.00000 1.0... m2 MULTILINESTRING ((0.50000 0.00000, 1.00000 2.0... dtype: geometry # 展示单个元素 ml[0]Polygon 面 在Polygon元素,其创建函数为: # shell表示单个多边形坐标,hole表示是否在内部创建空洞 Polygon(shell=None, holes=None) # 创建单个多边形 from shapely.geometry import Polygon pol = gpd.GeoSeries([Polygon([(0, 0), (4, 1), (5, 3), (0, 7)]), Polygon([(0, 0), (1, 1), (1, 2)])], index=['pol1','pol2']) pol pol1 POLYGON ((0.00000 0.00000, 4.00000 1.00000, 5.... pol2 POLYGON ((0.00000 0.00000, 1.00000 1.00000, 1.... dtype: geometry # 展示单个元素 pol[0]MultiPolygon 面集 # MultiPolygon([Polygon,Polygon])用于创建多个多边形的集合 from shapely.geometry import Polygon, MultiPolygon mpol = gpd.GeoSeries([MultiPolygon([Polygon([(0, 0), (0, 1), (1, 1), (1, 0)]), Polygon([(2, 3), (3, 2), (3, 3)])])]) mpol 0 MULTIPOLYGON (((0.00000 0.00000, 0.00000 1.000... dtype: geometry # 展示单个元素 mpol[0]GeoSeries常用属性如下: area: 计算GeoSeries每一个元素的面积,面积单位取决于坐标系 bounds: 获得GeoSeries每一个元素的最小外接矩形的左上角坐标和右下角坐标 total_bounds:计算整个GeoSeries对象的最小外接矩形的左上角坐标和右下角坐标 length:计算GeoSeries每一个元素的周长,周长单位取决于坐标系 geom_type: 获得GeoSeries每一个元素的类型 is_valid: 判断GeoSeries每一个元素是否构成合理的多边形 centroid:获得GeoSeries每一个元素的几何中心 # 创建包含两个多边形的GeoSeries from shapely.geometry import Polygon pol = gpd.GeoSeries([Polygon([(0, 0), (4, 1), (5, 3), (0, 7)]), Polygon([(0, 0), (1, 1), (1, 8)])], index=['pol1','pol2']) # 计算各个元素的面积 pol.area pol1 21.0 pol2 3.5 dtype: float64 # 计算各元素的最小外接矩形的顶点数据 pol.bounds minx miny maxx maxy pol1 0.0 0.0 5.0 7.0 pol2 0.0 0.0 1.0 8.0 # 计算整个对象的最小外接矩形的顶点数据 pol.total_bounds array([0., 0., 5., 8.]) # 计算GeoSeries每一个元素的周长 pol.length pol1 19.762298 pol2 16.476471 dtype: float64 # 获得GeoSeries每一个元素的类型 pol.type pol1 Polygon pol2 Polygon dtype: object # 判断GeoSeries每一个元素是否构成合理的多边形 pol.is_valid pol1 True pol2 True dtype: bool # 获得GeoSeries每一个元素的几何中心 pol.centroid pol1 POINT (1.88889 3.00000) pol2 POINT (0.66667 3.00000) dtype: geometryGeoSeries常用方法如下: distance(): 计算GeoSeries每一个元素到某一几何体的距离,具体使用查看官方文档 representative_point():返回位于GeoSeries每一个元素中的坐标点,该点不会是几何中心。 to_crs():改变GeoSeries的坐标参考系 plot():绘制整个GeoSeries # 创建包含两个多边形的GeoSeries from shapely.geometry import Polygon pol = gpd.GeoSeries([Polygon([(0, 0), (4, 1), (5, 3), (0, 7)]), Polygon([(0, 0), (1, 1), (1, 8)])], index=['pol1','pol2']) # 计算距离 point = Point(1, 0) pol.distance(point) pol1 0.242536 pol2 0.707107 dtype: float64 # 返回位于GeoSeries每一个元素中的坐标点 pol.representative_point() pol1 POINT (1.25000 5.00000) pol2 POINT (0.78125 4.50000) dtype: geometry # 绘制整个GeoSeries pol.plot() 1.2 GeoDataFrameGeoDataFrame是包含GeoSeries的表格数据结构,GeoDataFrame最重要的属性是它始终具有一个表示地理状态信息的GeoSeries列,即geometry列。调用GeoDataFrame的方法或属性时,将始终作用于geometry列。 下面用一个世界地图的实例介绍GeoDataFrame常用函数的使用方法,具体函数的接口请查阅官方文档。 ⚠⚠⚠注意该世界地图数据集来自geopandas自带数据集,部分数据数据可能缺失或错误,且数据集仅可用于学习不做其他用途 import geopandas as gpd # 加载自带世界地图数据集 world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres')) # 查看数据集 # 人口,大洲,国名,国家缩写,ISO国家代码,gdp,地理位置数据 world.head() pop_est continent name iso_a3 gdp_md_est geometry 0 920938 Oceania Fiji FJI 8374.0 MULTIPOLYGON (((180.00000 -16.06713, 180.00000... 1 53950935 Africa Tanzania TZA 150600.0 POLYGON ((33.90371 -0.95000, 34.07262 -1.05982... 2 603253 Africa W. Sahara ESH 906.5 POLYGON ((-8.66559 27.65643, -8.66512 27.58948... 3 35623680 North America Canada CAN 1674000.0 MULTIPOLYGON (((-122.84000 49.00000, -122.9742... 4 326625791 North America United States of America USA 18560000.0 MULTIPOLYGON (((-122.84000 49.00000, -120.0000... # 绘制GeoDataFrame # column选择为地图着色的数据列 # kind设置图片类型,某些类型会导致column参数不可用 # cmap设置matplotlib调色板 # legend设置是否显示图例 # figsize设置图片大小 # legend_kwds设置matplotlib的legend参数,该参数需要根据实际图例形式设置 world.plot(column="continent",kind="geo", cmap="coolwarm",legend=True,figsize=(9,6), legend_kwds={'loc': (0,1.02),'ncol':4})此外,pandas的数据操作也可以用于geopandas。 # 查看数据 world.iloc[:4,[3,4]] iso_a3 gdp_md_est 0 FJI 8374.0 1 TZA 150600.0 2 ESH 906.5 3 CAN 1674000.0 # 查看人口小于10万的国家地区 world.loc[world['pop_est']选择需要的数据格式->选择是否包含子区域->复制JSON API->read_file()函数读取文件即可。例如本文选择江苏省的地理信息,包含子区域的意思为是否包含江苏省内地区的地理数据。GeoPandas提供了各种数据集筛选函数,以提取感兴趣区域。其在0.1.0版本中新增了Bounding Box Filter,在0.7.0版本中新增了Geometry Filter和Row Filter,以及其他需要Fiona库的过滤器。 2.2.1 Bounding Box FilterGeoPandas在0.1.0版本引入了Bounding Box Filter。Bounding Box Filter用于提取与边界框相交的子区域,注意边界框格式为(左下角x, 左下角y, 右上角x, 右上角y)。 # 框选数据 import geopandas as gpd # 设置框选范围 bbox = (119, 32, 121, 33) # 读取数据 data = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/320000_full.json", bbox=bbox) # 结果是筛选与bbox相交的子区域 data.plot()GeoPandas在0.7.0版本引入了Geometry Filter。与Bounding Box Filter类似,Geometry Filter也是筛选与指定区域相交的子区域,但Geometry Filter输入的是几何图形。 # 查看几何图形与哪些子区域相交 from shapely import geometry # 设置框选范围 p = geometry.Polygon([(119,32),(120,33), (121, 33)]) # 读取数据 data = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/320000_full.json",mask=p) # 筛选区域 data.plot()GeoPandas在0.7.0版本引入了Row Filter。Row Filter用于提取指定行范围的数据。 # 读取前3行数据 data = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/320000_full.json",rows=3) data adcode name childrenNum level parent subFeatureIndex geometry 0 320100 南京市 11 city {'adcode': 320000} 0 MULTIPOLYGON (((119.1 32.5, 119.1 32.5, 119.1 ... 1 320200 无锡市 7 city {'adcode': 320000} 1 MULTIPOLYGON (((119.5 31.2, 119.5 31.1, 119.6 ... 2 320300 徐州市 10 city {'adcode': 320000} 2 MULTIPOLYGON (((118.4 34.4, 118.4 34.4, 118.4 ... # 读取第1行和第2行数据 data = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/320000_full.json",rows=slice(1, 3)) data adcode name childrenNum level parent subFeatureIndex geometry 0 320200 无锡市 7 city {'adcode': 320000} 1 MULTIPOLYGON (((119.5 31.2, 119.5 31.1, 119.6 ... 1 320300 徐州市 10 city {'adcode': 320000} 2 MULTIPOLYGON (((118.4 34.4, 118.4 34.4, 118.4 ... 2.2.4 cx索引器除了标准的Pandas方法外,GeoPandas还提供使用cx索引器进行基于坐标的索引。通过cx索引器可以挑选特定区域。 # 读取世界地图 world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres')) # 筛选与经度范围为80到200,纬度小于0相交的区域,这里坐标代表经纬度。 southern_world = world.cx[80:200, :0] southern_world.plot(figsize=(9, 6))地信、计算机和遥感等领域的从业人员或多或少都会接触地理信息系统(GIS,Geographic Information System)的相关知识。所谓GIS简单来说就是一个以计算机为核心,对地理空间位置相关数据进行创建、管理、分析、绘制和展示的多功能集成信息系统。绘制地图,则需要了解GIS中的坐标参考系。本文只是简单介绍坐标参考系的相关内容。相关内容总结于以下文章,想要具体了解更多内容也可以看看这些文章。 Coordinate Reference System and Spatial Projection 聊聊GIS中的坐标系 基于geopandas的空间数据分析——坐标参考系篇 前后端都要懂一点的 GIS GIS基础知识如下图所示,对于二维平面上的任何一点,我们可以通过如下的笛卡尔坐标系确定其唯一坐标,并用\((x,y)\)表示。 然而我们生活在一个恰好是“圆形”的三维地球上,一般来说要定义地球上物体的位置,需要一个适应地球形状的三维坐标系。但当我们在纸上或计算机屏幕上展示地球时,受限于存储介质(纸张、沙盘、书本、屏幕)原因,不得不使用二维平面来表达三维的地球。这种将三维球面转换为二维平面所用到的坐标系参考就是坐标参考系统(CRS,Coordinate Reference System),如下图所示。 常见情况下,根据不同的转换方式,CRS又可分为地理坐标系和投影坐标系。 地理坐标系统 参考平面为球面,用球面坐标来表示地球上的位置,叫做地理坐标系统。这些坐标系的坐标单位通常为角度,如经纬度。所谓经度就是指离被称为本初子午线的南北方向走线以东或以西的度数。纬度则是指地球某点与球心的连线和地球赤道面所成的线面角。如下图所示: 在数据表示中,东经正数,西经为负数,经度范围为-180度到180度。赤道以北为北纬,以南为南纬,纬度范围为-90度到90度。按照常用地理坐标系统表示的二维地图如下所示: 投影坐标系统 由于我们所在的地球实际不是规则球体,一经度和一纬度代表的实际距离是不相等的(赤道除外)。这导致我们无法通过经纬度坐标精确计算某块地区的面积,长宽。也无法精准表达某一块地区的实际形状。因此,有必要采取某些方法将地球表面投影到一个二维平面,这种地理投影方法也就是投影坐标系,实际作用为把地球表面经纬度按照一定的数学规则转换为平面坐标。 投影坐标系在二维平面中定义,其始终基于地理坐标系,位置由二维坐标标识,度量单位用米、千米等来表示。不同的投影方法会导致不同的投影结果,因此各种各样的投影坐标系统被提出。但是不管是哪种地图投影方法,都会导致地图投影变形。三维投影二维出现投影变形是不可避免的,但不会在任何位置都会产生变形。那些没有发生变形的点或线称为标准点线,距离标准点线越近,变形程度则越小。因此在实际投影某块区域时,最好选择标标准点线离该区域近的投影坐标系。 关于常见投影方法的介绍可以阅读本节开头给出的参考文章。想要感受不同投影方法带来的效果,可以使用一个非常酷的坐标系投影介绍网站:worldmapcreator。在worldmapcreator中,我们可以选择不同的投影方法来实时显示投影变换效果,以及自定义投影结果,具体如下: EPSG:European Petroleum Survey Group(欧洲石油调查小组)是一个涉及测地学、测量、制图学与石油勘探相关的科学组织,它维护和发布的EPSG编码系统提供了各种坐标参照系统crs的数据集参数。我们通过如下网站查询EPSG投影坐标的编码:espg.io,该网站界面介绍如下: 互联网常见的坐标参考系如下图所示。WGS84和CGCS2000都为地理坐标系统。其中WGS84(World Geodetic System 1984)是目前应用最为广泛的坐标系统。在没有特别声明的时候,默认使用WGS84坐标系统存储地图数据或提供地图数据。WGS84坐标单位为度,WGS84的EPSG代码为EPSG:4326。CGCS2000是我国当前最新的国家大地坐标系,它的EPSG编码为EPSG:4490,在大多数精度要求不高的情况下,我们可以默认WGS84坐标系等同于CGCS2000坐标系。 除此之外,WGS84 Web Mercator(EPSG:3857)坐标系统也较为常用。EPSG:3857为投影坐标系统,单位为米。EPSG:3857在WGS84坐标系基础上进行伪墨卡托投影(Pseudo-Mercator)得到。关于Pseudo-Mercator介绍见墨卡托投影。墨卡托投影能够保证方向的正确,所以EPSG:3857广泛用于导航海航。但越到高纬度,墨卡托投影会导致大小扭曲越严重,到两极会被放到无限大。在长度表示上EPSG:3857有很大偏差,无法用于实际距离和面积的测算。墨卡托投影下每个国家的大小和实际大小的差异见下图: 当然CRS还有其他常见的数据编码格式。计算机领域更多地使用EPSG编码。 3.3 在GeoPandas中管理坐标参考系GeoPandas提供了丰富的函数接口以管理坐标参考系,具体使用如下: 3.3.1 坐标参考系管理查看坐标管理系 # 框选数据 import geopandas as gpd # 读取江苏省数据 data = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/320000_full.json") # 查看坐标管理系,可以看到默认坐标参考系为EPSG:4326,单位为度。 data.crs Name: WGS 84 Axis Info [ellipsoidal]: - Lat[north]: Geodetic latitude (degree) - Lon[east]: Geodetic longitude (degree) Area of Use: - name: World. - bounds: (-180.0, -90.0, 180.0, 90.0) Datum: World Geodetic System 1984 ensemble - Ellipsoid: WGS 84 - Prime Meridian: Greenwich # 计算面积,可以看到面积计算结果不对,因为单位是度 # 直接计算面积会报错,提示UserWarning: Geometry is in a geographic CRS sum(data.area) 9.997823062190006 # 绘图 data.plot()设置坐标参考系 对于有些几何形状或者坐标点没有设置坐标参考系,为了告诉GeoPandas如何解释这些几何形状,需要设置默认坐标系。GeoPandas提供set_crs函数来执行这一操作。 p = gpd.GeoSeries([geometry.Point([134.1233,23.456])]) # 没有默认坐标系 p.crs # 重新设置坐标系 p = p.set_crs("EPSG:4326") p.crs Name: WGS 84 Axis Info [ellipsoidal]: - Lat[north]: Geodetic latitude (degree) - Lon[east]: Geodetic longitude (degree) Area of Use: - name: World. - bounds: (-180.0, -90.0, 180.0, 90.0) Datum: World Geodetic System 1984 ensemble - Ellipsoid: WGS 84 - Prime Meridian: Greenwich重设坐标参考系 GeoPandas提供了crs函数以重新设置坐标参考系。 world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres')) world.crs Name: WGS 84 Axis Info [ellipsoidal]: - Lat[north]: Geodetic latitude (degree) - Lon[east]: Geodetic longitude (degree) Area of Use: - name: World. - bounds: (-180.0, -90.0, 180.0, 90.0) Datum: World Geodetic System 1984 ensemble - Ellipsoid: WGS 84 - Prime Meridian: Greenwich world.plot()坐标参考系选择方法 如果我们获得是地理坐标系数据,如WGS84,在没有更多信息的情况下,如何转换到投影坐标参考系。一个好的办法就是根据所在区域的几何中心点,从前面提到的espg.io寻找合适EPSG编码。这种方法很不专业,只是低精度使用! 以投影江苏省地图为例,一种办法是找到整个国家的几何中心,一种是找到当地的几何中心。 # 框选数据 import geopandas as gpd # 读取整个国家的数据,无子区域数据 data_china = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/100000.json") # 查看经纬度几何中心,POINT (103.9 36.4) # 东经103.9 data_china.centroid 0 POINT (103.9 36.4) dtype: geometry # 框选数据 import geopandas as gpd # 读取江苏省数据,无子区域数据 data = gpd.read_file("https://geo.datav.aliyun.com/areas_v3/bound/320000.json") # 查看经纬度几何中心,POINT (119.5 33.0) # 东经119.5 data.centroid 0 POINT (119.5 33.0) dtype: geometry接下来我们可以在 espg.io中,搜索New Beijing kind:PROJCRS。这个搜索词表示搜索以New Beijing为关键词,投影坐标(kind:PROJCRS)的结果。当然以beijing kind:PROJCRS和xian kind:PROJCRS也是可以的。这些都是我国常用的坐标系,不添加kind:PROJCRS也没关系,搜索结果Coordinate reference systems中有相应的选项。然后我们查看地区的几何中心在是否坐标系简介的Area of use中。如下图所示: 比如以整个中国的几何中心(103.9 36.4)查找的结果为EPSG:4573,如下所示: 以江苏省的几何中心(119.5 33.0)查找的结果为EPSG:4586,如下所示: 坐标参考系转换实例 江苏省的面积为10.72万平方公里,通过EPSG:4573计算的面积为10.84万平方公里,EPSG:4586计算的面积为10.38万平方公里。这两种都是近似估算,实际不会这么粗糙。但是可以看到选择合适的EPSG编码结果都大差不差,想要更准确的结果,可以自行查找其他EPSG编码。 此外,为了展示结果可靠性 我们可以获得某一坐标位置进行绘制。通过地图查询苏州市上的某一经纬度坐标(120.601868 31.332525),如下所示: 从上面可以看到EPSG:4573和EPSG:4586表示地图的横坐标单位不一样。EPSG:4573代表的是New Beijing / Gauss-Kruger zone 18,这个zone就是指带号的意思。EPSG:4586代表的是New Beijing / Gauss-Kruger CM 117E,表示没有加带号。简单解释:例如常见的6度分带投影,表示从零度子午线开始,自西向东每个经差6度为一投影带,全球共分60个带,用1,2,3,4,5,... ,60表示。比如没带号的横坐标为500000米。有带号zone 18,表示在没带号的横坐标前加上带号18,变为18500000。具体可以查阅投影分带。 在本例中,如果以江苏省的几何中心查找加带号的坐标系,则应该是New Beijing / Gauss-Kruger zone 20(EPSG:4575)。与EPSG:4586相比,以EPSG:4575为基准绘图,横坐标加了带号(20),面积计算结果和是相等的。这些只是GIS基础知识,本文说得非常简单也有诸多错误(因为作者非专业人士),低精度软件绘制地图也不需要用到过多专业的知识。如果想高精度绘图,请从头开始学习GIS相关知识。 data.to_crs("EPSG:4575").area.sum()/1e6/1e4 10.381446380451155 ax = data.to_crs("EPSG:4575").plot() p = gpd.GeoSeries([geometry.Point([120.601868,31.332525])]).set_crs("EPSG:4326") p.to_crs("EPSG:4575").plot(ax=ax, color='red', markersize=100) # EPSG:4586 为POINT (842904.5 3473513.1) p.to_crs("EPSG:4575") 0 POINT (20842904.5 3473513.1) dtype: geometry本文来自博客园,作者:落痕的寒假,转载请注明原文链接:https://www.cnblogs.com/luohenyueji/p/17222011.html ![]() |
今日新闻 |
推荐新闻 |
CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3 |